################################################################
#
#  Comparação da geração de normais com Box-Muller e gerador dos R
#
###############################################################

boxmuller <- function(){
 n <- 10^4
 samples <- matrix(ncol=2,nrow=n)
 set.seed(1)
 for (i in 1:n){
 u1 <- runif(1)
 u2 <- runif(1)
 R <- sqrt(-2*log(u1))
 theta <- 2*pi*u2
z1 <- R*cos(theta)
z2 <- R*sin(theta)
 samples[i,1] <- z1
 samples[i,2] <- z2
 }
 norm01 <- rnorm(n)
 geração <- rep(c("z1","z2","normal"),n)
 value <- c(samples[,1],samples[,2],norm01)
 df <- data.frame(value,geração)
 library(ggplot2)
 plt <- ggplot(df, aes(x=value, color=geração, fill=geração)) +
 geom_histogram(aes(y=..density..), bins = 60, position="identity", alpha=0.3) + labs(x="Value", y="Density") +
 theme_bw()
 print(plt)
 }
 boxmuller()

######################################################################################################

set.seed(1)
 n <- 10^4
 samples <- matrix(ncol=2,nrow=n)
 for (i in 1:n){
 u1 <- runif(1)
 u2 <- runif(1)
 R <- sqrt(-2*log(u1))
 theta <- 2*pi*u2
z1 <- R*cos(theta)
z2 <- R*sin(theta)
 samples[i,1] <- z1
 samples[i,2] <- z2
 }

 norm01 <- rnorm(n)
 generation1<- rep(c("z1","normal"),n)
 value1 <- c(samples[,1],norm01)
 generation2<- rep(c("z2","normal"),n)
 value2 <- c(samples[,2],norm01)
compare_data1 <- data.frame(value1,generation1)
compare_data2 <- data.frame(value2,generation2)
dx <- seq(-4, 4, by = 0.01)
curves <- data.frame(x = dnorm(dx)) 

library(ggplot2)
library(patchwork)

 


 p1 <- ggplot(compare_data1)+
 geom_histogram(aes(x=value1, y = stat(density)), fill = "blue", binwidth = 0.1, bins = 60, position="identity", alpha=0.3) + 
stat_function(fun = dnorm, args = list(mean = 0, sd = 1))+
labs(x="generation1", y="Density") +
 theme_bw()

 p2 <- ggplot(compare_data2)+
 geom_histogram(aes(x=value2, y = stat(density)), fill = "green", binwidth = 0.1, bins = 60, position="identity", alpha=0.3) + 
stat_function(fun = dnorm, args = list(mean = 0, sd = 1))+
labs(x="generation2", y="Density") +
 theme_bw()
p1 + p2


#####################################################################################
#
#  Comparação da geração de normais com: Transformação de Box-Muller vs Soma de 12 U(0,1) indep.
######################################################################################

nsim=10000
u1=runif(nsim)
u2=runif(nsim)
X1=sqrt(-2*log(u1))*cos(2*pi*u2)
X2=sqrt(-2*log(u1))*sin(2*pi*u2)
U=array(0,dim=c(nsim,1))
for(i in 1:nsim)U[i]=sum(runif(12,-.5,.5))
par(mfrow=c(1,2))
hist(X1)
hist(U)
a=3
mean(X1>a)
mean(U>a)

compare_data <- data.frame(X1,U)

library(ggplot2)
library(patchwork)

p1 <- ggplot(data=compare_data) + geom_histogram(aes(x=X1),binwidth=0.2)
p2 <- ggplot(data=compare_data) + geom_histogram(aes(x=U),binwidth=0.2)

p1 + p2

#######################################################################